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Abstract 

The growth in number and nature of dynamical at- 
tractors in Kauffman NK network models are still not 
well understood properties of these important random 
boolean networks. Structural circuits in the underpin- 
ning graph give insights into the number and length dis- 
tribution of attractors in the NK model. We use a fast 
direct circuit enumeration algorithm to study the NK 
model and determine the growth behaviour of struc- 
tural circuits. This leads to an explanation and lower 
bound on the growth properties and the number of at- 
tractor loops and a possible K-relationship for circuit 
number growth with network size N. We also in- 
troduce a mixed-K model that allows us to explore 
N (K) between pairs of integer K values in Kauffman- 
like systems. We find that the circuits' behaviour is 
a useful metric in identifying phase transitional be- 
haviour around the critical connectivity in that model 
too. We identify an intermediate phase transition in cir- 
cuit growth behaviour at K = Ks ~ 1.5, that is dis- 
tinct from both the percolation transition at Kp = 1 
and the Kauffman transition at Kq = 2. We relate this 
transition to mutual node reachability within the giant 
component of nodes. 



1 Introduction 



as quantum gravity through their relationship with </> 3 - 
networks [5, 6]. RBNs have many interesting proper- 
ties [7] and have been amenable to various analyses [8] 
including mean-field theory. They also continue to be 
an important and interesting tool in studying biological 
and artificial life problems [9, 10]. 



NK Kauffman Network Model Phase Diagram 
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Figure 1 : Phase diagram for the Kauffman NK Model 
in terms of (integer- valued) connectivity K and bit-bias 

p. 



Kauffman 's NK-Model [1,2] of an TV-node random 
boolean network with if -inputs to a boolean function 
residing on each node has found a significant role in the 
study of complex network properties. Random Boolean 
Network (RBN) models are effectively a generalisation 
of the 1 -dimensional Cellular Automata model [3] and 
have important applications in biological gene regula- 
tory networks [4] but also in more diverse areas such 



One key property of RBNs is the now well estab- 
lished existence of a frozen phase and a chaotic phase 
[11, 12] and the critical phase transition lies at the in- 
teger value of connectivity K c = 2 ,*_ , = 2 for un- 
biased networks with a mean boolean function output 
value of p = 0.5. This gives rise to the phase diagram 
shown in figure Q] which is dominated by the chaotic 



region and for which simulations are limited to the in- 
teger K values shown. It has therefore been of most 
interest to study RBNs at or around this critical (inte- 
ger) value of K c = 2. In this paper we explore another 
mechanism to explore the critical region by adopting a 
mixed- K system. 

The Random Boolean Network or graph G is ex- 
pressed as a four-tuple G = (V, E, F, x) and has N = 
\V\ = \F\ = \x\ nodes or vertices^), and N a = \E\ 
directed edges or arcs, which express the Ki inputs for 
node i. Formally we let F2 = {0, 1} be a bit-field or 
Galois Field with the operator defined modulo 2 and 
¥2 is the set of all possible vectors of K bits. Vector 
x G ¥2 and X{ is then the bit at position i. 

The Kauffman NK-Model Network is constructed 
with fixed (integer) K = 1,2,3,.. and a boolean func- 
tion fi G F of Ki inputs is assigned to each node. All 
the nodes of the network carry a boolean variable xi 
which may be initialised randomly and each of which is 
updated (usually, but not necessarily) synchronously in 
time t from its j -labelled inputs so that: 

Xi(t) ^ fi ( Xj (t -1)) J = 1,2,..., Ki (1) 

The boolean functions / thus map Ff" 1— > F2. The map- 
ping can be expressed as a truth table and readily imple- 
mented as look-up table in a simulation program [13]. 

Work has been carried out on a number of different 
time-update mechanisms for boolean networks includ- 
ing asynchronous algorithms [14, 15]. In this paper we 
consider only synchronous updates where all nodes ex- 
ecute their boolean function once, at the same time, at 
every time step. Other studies have also considered how 
noise [16] effect the crossing times between distinct tra- 
jectories. In this paper we adhere to the quenched con- 
vention whereby node connections are assigned once 
and for all time, and particular boolean functions are 
assigned randomly and uniformly to nodes once and for 
all time. 

The NK-network model assigns the Ki inputs for 
node i randomly and with uniform probability across 
all possible nodes. Even for a large network there is 
still a non-zero probability of assigning a node as one 
of its own inputs. In the case of Ki > 1 there is also a 
possibility of assigning a node j as an input of i more 
than once. These self-edges and multiple-edges can 
have a subtle effect on the behaviour of the NK-network 
model [17]. 

A significant body of work has now been carried out 
on the roles of different sub-classes of boolean func- 
tions including the so called canalizing functions [18] 
and in particular the effect of bit-bias thresholds and 



frozen or fixed- value boolean functions on particular el- 
ements of the network [19]. A network can therefore 
be restricted to only have some subset of the possible 
boolean functions. In the work we describe here, we 
use network sizes large enough to sample all possible 
if -input boolean functions for the largest K present in 
the system. 

An important consequence of the boolean functions 
in RBNs is the formation of attractor basins [20] . These 
are observed in RBN models whereby diverse initial 
starting conditions will still lead to statistically similar 
behaviour. The state of the network falls into attrac- 
tor cycles whereby a chain of interdependence of nodes 
(via their boolean functions) leads to the network pe- 
riodically repeating its state. The number and length 
of these periods or attractors is of great importance in 
understanding the behaviour of the NK-model and as- 
sociated application problems. This can be seen quan- 
titatively by tracking a metric such as changes in the 
normalised Hamming distance between the network's 
successive boolean states. We discuss this metric in sec- 
tional 

Of particular recent interest in the literature has been 
the uncertainty concerning the number of attractors [2 1 , 
22] and how their number and lengths varied with the 
size of the network. Scaling was initially believed to be 
0{y/~N) [23]. It was later reported as linear [24], and 
then as "faster than linear" [25] and subsequently as 
"stretched exponential" in [26, 27] but is now known to 
be faster than any power law [28]. 

A recent review of the RBN model [8] discusses 
the attractor behaviour in terms of the loops of boolean 
variable states that form and several exact results con- 
cerning these loops have been obtained for the case of 
connectivity K = 1 [29] . Important observations con- 
cern the distribution of components with particular sub- 
classes of possible boolean functions. These "relevant 
elements" are defined as those nodes that are not frozen 
and that control at least one other relevant element in 
the system [27]. A number of important results have 
been obtained using particular sub-classes of the pos- 
sible boolean functions. Drossel et al. [22] have con- 
sidered networks with non-fixed boolean functions thus 
making all elements relevant and have therefore shown 
the equivalence of K = 1 and K = 2 networks under 
appropriate restrictions on the boolean functions. 

In this paper we use numerical methods to investi- 
gate the role that structural circuits play in the complex 
structure of the network and the resulting attractor be- 
haviour of RBNs both for mono-if and mixed-if net- 
work systems. 




Figure 2: 12 Node Network, with (K) w 2.5. For each 
node, input degree is shown on the left and output de- 
gree is shown on the right. One particular circuit, 1-10- 
7-8-1, is shown. 



Recent work in the literature has used trajectory sam- 
pling to study attractor behaviour. The combinatorics of 
RBN models means that the number of boolean func- 
tions grows as 2 2 with a consequent rapid growth in 
the number of possible network states with network 
size. Taking limited numbers of sample trajectories 
through this state space can lead to very misleading 
results. Numerical sub- sampling of attractor trajecto- 
ries seems to be the main difficulty behind obtaining a 
good understanding of attractor scaling. In this paper 
we explore the structural properties of RBNs including 
the number and length distribution of elementary cir- 
cuits and of components. We compute these properties 
exactly using brute force enumeration techniques for a 
range of network sizes and connectivities. Our statisti- 
cal sampling is only over different randomly configured 
networks, not over attractor trajectories. 

We have found it necessary to study quite large 
samples (up to 100,000) of networks of size up to 
250,000 nodes. The network size N must be at least 
large enough to adequately sample all possible if-input 
boolean functions for the largest K present in our mixed 
systems. There are some good software tools avail- 
able for experimenting with RBNs such as those of Ger- 
shenson [7] and Wuensche's DDLab [30], but we were 
forced to develop our own custom D code to simulate 
very large-scale systems [13]. 

The so-called "edge of chaos" regime [2] is quite 
narrow when K is restricted to integer values. In this 
paper we also consider other ways to explore the model 
phase space. We explore a mixed-if system which we 
refer to as the N (K) Network model, with the under- 
standing that although individual nodes (must) have an 
integer K number of inputs, the system as a whole can 
have a mean, or effective, K value if there is a distribu- 
tion of nodes each with different number K of inputs. 
It then becomes a matter of deciding on a sensible Re- 
distribution for a given model system. 

Although some work has been done developing sim- 
ulations that employ mixed-K models with Poisson or 
other distributions [30] we believe no one has yet stud- 
ied mixed models in systematic detail and furthermore, 
that our pair-wise model is a novel way of achieving 
a definite i^ e ff. In [31] Skarja et al. employed a 
skewed binomial distribution of K values and found an 
enhanced tendency to orderliness (stability) but based 
their study of attractors on trajectory sampling. 

It is not altogether obvious what the meaning of a 
distribution in K might mean if it has a long tail with 
values at high and low K relative to the mean, such 
as the Poisson distribution in node outputs that results 



from the mono-K NK model. Consequently we have 
chosen to study a pair-wise model that generates a K 
value for each node that interpolates between a pair of 
integer K values. We define a parameter P that is the 
linear probability of a node having K = K 2 rather than 
K — K\. So the cases P = or P = 1 correspond 
to the pure integer states of all nodes having K 2 or K\ 
respectively, and P = 0.5 constitutes an equal mix of 
the designated pair of K values. This means we can ap- 
proach the transitional value of K = 2 from above or 
below and also means we can explore properties with In- 
dependent relationships thoroughly over a larger range 
of fractional K values rather than just the small subset 
of (small) integer K values that is feasible numerically. 

Figure [2] shows a small network generated with our 
mixed pair model. Nodes have 2 or 3 inputs but can 
have zero or more outputs depending upon the random 
distribution used during construction in the normal way 
for NK networks. 

In [17] we explored the behaviour of the growth of 
the number of elementary circuits and the resulting cir- 
cuit length distributions in integer if -valued NK net- 
works. We suggested that the number of circuits gave a 
lower bound on the number of attractors present in the 
associated NK network model. Our numerical exper- 
iments exactly enumerated the circuits in various NK 
networks, without resorting to trajectory sampling and 
we concluded that the growth in the number of attrac- 
tors had to be at least as fast as exponential in network 
size N. This clarified some of the recent controversy 
concerning the growth rates in the number of attractors. 

In our earlier work we had an insufficient range of K 
values to determine a circuit growth relationship with 
K. In this present paper we have sampled higher K 
values and also examined the pair- wise or fractional K 
model over a large range of K values and have therefore 
been able to determine likely relationships between K 
and the growth of the number of circuits with network 
size N. 

We have also observed long-time variations in the 
normalised Hamming distance for various K-valued 
networks at small and large network sizes. These sug- 
gest the very strong importance of adequately sampling 
all possible Boolean functions for a given K value. This 
phenomenon may also explain the anomalous and mis- 
leading results on attractor growth obtained from sam- 
pling trajectories in too- small networks. 

In section [2] we discuss the pair- wise K model and 
implementation issues. We summarise the the role cir- 
cuits appear to play in Kauffman nets in section [3] and 
their enumeration in section [4] In section [5] we present 



results on the number of circuits and their length distri- 
bution both for integer and fractional K network sys- 
tems. In sections [6] and [7] we offer some discussion of 
the results and conclusions concerning circuit growth 
with K and TV and the properties of the mixed-pair K 
model. 



2 K Pairs and the TV (K) Model 

The conventional Kauffman NK model with mono-K 
can be extended to a system with a distribution of K- 
input nodes. RBN simulation tools like Wuensche's 
DDLab [30] do make provision for this but it is unclear 
how to systematically investigate possible distributions, 
particularly when large samples are required. We might 
intuit that the mean or effective K value for the whole 
system plays an important role, but it is not clear how 
smeared-out any behaviour might be that results from 
multiple K values. Our pair- wise model approaches the 
critical value in K from either side by adopting a sim- 
ple uniform mix of just two possible values of K. Most 
useful are if = 2, 3 and K = 1, 2 to approach K c = 2 
from above or below, respectively. 

An effective-if for the whole network can be defined 
as: 

*eff = (Ki) = (1 " P)-K low + RK high (2) 

where any individual node i has an integer valued num- 
ber of exactly K \ inputs and the average is over all N 
nodes. Individual nodes are randomly assigned (once 
and for all time) their particular K value. 

We can investigate both static and dynamic proper- 
ties of this simple mixed model. Static properties are 
measured from simple graph analysis, and dynamics 
can be obtained by examining the Normalised Ham- 
ming distance between subsequent bit states of the net- 
works' nodes. 

If we have the vector x € Ff" and X{ is the bit at po- 
sition i in the network, we define the Hamming weight 
as 

w H = \{i\xi^0,i=l,2,...,K}\ (3) 

and the Hamming distance between two vectors x, y as 



H d (x,y) =w H (x®y) 



(4) 



It is useful to normalise this by the network size N. 
Of particular interest is the Hamming distance between 
subsequent states of the network. This is easily calcu- 
lated as 1 — asame where asame is the fraction of nodes 



1-Step Lag Correlation in Kauffman Networks (256,000 Nodes) 
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Figure 3: Single-Step lag correlation function (or frac- 
tion asame) for 256, 000 Node Kauffman Networks for 
if = 0,l,2,3,4,5 and networks with mixed pairs of K 
values in between. This is the same as "one minus the 
normalised Hamming distance." 



that have the same bit value at subsequent steps or the 
single- step correlation function. 

Figure[3] shows the single step correlation function as 
measured from a variety of different K values. At low- 
K the network is barely connected and it very rapidly 
converges to a fixed Boolean state. In highly connected 
networks the attractor loops (of various lengths) intro- 
duce periodic cycles of correspondingly varied frequen- 
cies. Nevertheless the mean state of the network still 
converges to a stable value that depends critically upon 
K. 

Figure shows the long-term mean convergent val- 
ues of the Normalised Hamming distance lim^oo H\ 
between subsequent (ie time difference of 1) Boolean 
states of the network as it varies with mean K in our 
N (K) model. These are measured for 10 separate sam- 
ples of an N = 256, 000 node network. We observed 
that a stable reproducible long-term convergence value 
is obtained for each K- valued network only if N is sig- 
nificantly larger than the corresponding number Ns f = 
2 2 of possible Boolean functions. Not unreasonably, 
if TV is smaller than Ns f , then the network instantiation 
is not adequately sampling the Boolean function space. 
Empirically we found N > ANb s suffices to ensure a 
reasonable sample. Individual networks generally con- 
verge to a stable value within a few hundred time steps. 
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Figure 4: Normalised Hamming long term mean values 



To avoid possible lack of convergence we discarded the 
first N time steps and averaged over a subsequent N 
steps. The average has to be taken over a time larger 
than any likely attractor traversal times present. 

Figure |4] does show a very clear transition at K = 2 
as expected for the integer-K NK model, but also in- 
dicates that the normalised Hamming distance changes 
smoothly for the N (K) model as we vary K. This sug- 
gests our pair- wise model is indeed a useful way to in- 
vestigate meaningful non-integer K- valued systems. In- 
tuitively we might expect that H\ (K — K c ) might have 
a straight-forward power or other closed function rela- 
tionship with (K), but we were not able to obtain a good 
numerical fit to any obvious forms. 



3 Structural Circuits and Attrac- 
tors 

Figure [2 shows a 12-node mixed-K network with K{ G 
2, 3, i^ e ff = 2.5 and where the construction algorithm 
has allowed self-arcs - in other words the inputs for each 
node have been chosen according to a flat uniform dis- 
tribution so they can connect to themselves. The conse- 
quent self-edges allow self-inputs in the corresponding 
RBN. These are known to play a vital role in supporting 
the number of attractors. A self-input or "self-ancestor" 
in the input dependence chain of boolean variables an- 
chors the periodic or attractor behaviour [8] of RBNs. 
We felt intuitively that the presence of structural cir- 



cuits would also be vital to the periodic attractor be- 
haviour and as Drossel et al. have shown there are defi- 
nite relationships between the number of attractors and 
the number of loops. In fact, the number of structural 
circuits provides a lower bound on the number of pos- 
sible attractors. Consideration of the exact number of 
enumerated circuits, and their distribution, gives insight 
into the controversy over the number of attractors in 
RBNs. 

An elementary circuit is a closed path along a subset 
of the edges of the graph such that no node, apart from 
the first and last, appears twice. The number of elemen- 
tary circuits for a fully connected graph is bounded by 



N-l 



C»_ i+1 (N-i)\ (5) 



as given by Harary [32]. This expression therefore rep- 
resents the limit for the number of structural circuits in 
an NK-network when K — > TV. 

Figure [2] shows one such circuit or loop in the net- 
work structure. In fact, exact enumeration (as shown in 
figure[5]) indicates that there are 30 arcs and 28 circuits, 
the longest of which is of length 10. It has 4 self- arcs 
(and hence two circuits are duplicated) and 3 multiple 
arcs. The maximum number of outputs is 5 and the min- 
imum is zero. If self-edges are disallowed we would ob- 
tain a higher number of circuits present in the network. 

4 Circuit Enumeration 

Various algorithms have been formulated to count the 
circuits in a graph but these either use infeasible 
amounts of memory or are time exponential [33, 34] 
with a time bound of 

0{N.e{c + 1)) (6) 

We count circuits using a variation of Johnson's algo- Fi S ure 5: Enumerated circuits for the network shown in 

rithm [35] implemented in D. For graphs of N vertices, figure |2J 

e edges, c circuits and 1 fully connected component, 
Johnson's algorithm is time bounded in time by 

0((N + e)(c+l)) (7) 

and space bounded by 0(N + e). Unlike Johnson's al- 
gorithm our code copes with partially connected graphs 
without resorting to the need to treat each of the possi- 
ble N c > 1 components separately [36]. This is still a 
highly expensive process since the number of circuits c 
itself grows very rapidly with (TV, e). 
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In the graph literature the term loop is unfortunately 
sometimes used to describe a self-edge or a circuit of 
length 1 . In the NK-networks we study the number of 
self-edges is much less than N, even for low K. How- 
ever we do count them and observe the effect of allow- 
ing them in the number of possible circuits and their 
length histograms. We have extended Johnson's pub- 
lished algorithm to cope with graphs with directed arcs 
(and not just bi-directional edges); with multiple com- 
ponents; and self- arcs. The computational complexity 
is not changed although we store some additional book- 
keeping information to support arcs. At worst this dou- 
bles the memory space required. All the work we re- 
port in this present paper is compute-time bound and 
not space bound. 

On a modern (circa 2007) compute server with 
4GBytes memory and a speed of 2.66GHz, we found 
it was entirely feasible to enumerate circuits exactly in 
networks of up to N « 140 for K = 1,2. Smaller 
networks were required for higher K. We were able to 
count components quite easily up to networks of around 
TV ^ 20, 000. We were able to exploit the near-linear 
speed-up of parallel job-farming to average our exact 
enumeration/counting results over many independently 
generated networks. 

In a detailed investigation of elementary circuits in 
the graph structure of NK networks [17] we found nu- 
merical evidence for rapid growth of the number of cir- 
cuits with network size N, but were unable to determine 
a reliable numerical relationship for growth in terms of 
K. This is largely because with integer K we are re- 
stricted to only a very few practical values. It is only 
feasible to study networks up to around K « 9. This 
situation is not likely to change even with linear im- 
provement in computer speeds or other supercomputing 
techniques since the growth in the number of circuits 
with K scales so rapidly. 

Using the mixed model however, we are able to in- 
vestigate intermediate values in K- space and attempt to 
find an empirical relationship for circuit growth with N 
for different mean values of K. 



Number of Circuits in N<K> Networks, 500 Samples 
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Figure 6: Number of circuits in Mixed-K Networks, 
sampled over 500 networks. 



5 Circuit Measurement Results 



In this section we present results of various numerical 
experiments, exactly enumerating the circuits over in- 
dependently generated sample networks. We emphasise 
that these data are not based on sampling attractor tra- 
jectories. 

Figure [6] shows the number of circuits in mixed-K 



networks, sampled over 500 different networks. Above 
a K value of 1.5 a straight line fit to \ogN c vs N is 
a good model for the data, whereas for low K, there 
appears to be a linear relationship between log N c and 
log TV, as shown in figure [71 
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Figure 7: Number of circuits in Mixed-K Networks, 
sampled over 500 networks. 



As we conjectured in [17] there is a definite change 
in behaviour between the low-K and high-K regimes, 
however this transition is not a simple one arising from 
the percolation transition at K = 1, as figure [7] indi- 
cates. 

We can explore this effect further by examining the 
least-squares fitted slopes obtained from figure [6] 

Figured shows a straight line fit to these fitted slopes 
from the data in figure [6] suggesting that the relationship 
N c « 2 s - N holds where s = s(K) and for K > 1.5 this 
is well fitted by a straight line so that s « 0.340i\~ and 
hence N c w A C 2°- M0K - N holds. 

Discounting the effect of self-arcs which grow in 
number linearly with N, and of multiple-arcs which 
grow with TV 2 , the mean number of connections Na is 
approximately K.N and therefore above K = 1.5 the 
model data strongly supports growth in the number of 
circuits with Na 

However, figure [8] indicates a very clear departure 
from this behaviour at K « 1.5. Interestingly the 
integer-K NK Kauffman model exhibits a transition in 
the long term Normalised Hamming distance at K — 2, 
as indeed does our mixed-K model as shown in figureH] 

As the error bars in figure [8] show the data does sup- 



port a straight line fit, although one could convince one- 
self there is a small anomaly at the critical K = K c = 
2. This is not surprising given the pairwise nature of 
our mixed model. We are essentially approaching K c 
independently from above and below. In the case of 
approach from above we have a system whose nodes 
mostly have K = 2 with a few of K = 3 whereas from 
below we mix in a minority with K — 1. It is in fact 
perhaps a point in favour of the simple pair- wise model 
that the two curves meet so closely. 



log_2 nCircuits vs N Slopes 
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Figure 8: Circuit Slopes of log2nCircuits vs N; fitted 
gradient is 0.340. Experimental data follows the curve 
shown, with a dark fitted straight line. Error bars were 
computed from standard deviations and are as shown. 



It is also instructive to examine the distribution of 
circuit lengths present in a typical system. 

Figure shows the length distribution of circuits in 
networks of fixed N but different K. At K values above 
1 .5, the distribution of circuit lengths has a definite peak 
around Af/2 and indicates that there is a non zero, how- 
ever small, possibility of Hamiltonian or near Hamilto- 
nian circuits present in the system that include all N 
nodes. Below K — 1.5 however the circuits length 
distribution falls monotonically with length and there 
is no modal length. Although we can only study rel- 
atively small network sizes in detail, we might expect 
that the fall off means there would be almost no cir- 
cuits of length greater than N/2 in large networks, and 
certainly, that the probability of there being any Hamil- 
tonian circuits is vanishingly small. 

Figure [TO] shows the modal circuit lengths as they 
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Figure 9: Distribution of circuit lengths in mixed K sys- 
tems, for fixed network size N = 40. 
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Figure 10: Modal Circuit lengths for 10,000 N = 40 
sample networks as it varies with i^ e ff . Note separate 
curves for networks with and without self- arcs. 



vary with effective K value. Above K = 1.5 the modal 
circuit length is non-trivial and grows logarithmically 
with mean K. Below the transition the most likely cir- 
cuit length present is unity for systems that allow self 
arcs, and two for systems that do not. 



6 Discussion on Network Compo- 
sition 

We have also carried out some standard graph metric 
analyses on our mixed system to clarify the role of the 
percolation transition at K = 1 on the growth circuits. 
Theoretically we expect that, ignoring the effect of self- 
arcs and multiple arcs, the percolation transition for in- 
finite sized networks is exactly K = 1 [37]. 

Number of Components in NK Networks, 1000 Samples 
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Figure 11: The number of cluster components in a 
mixed K system for different K and network sizes. 



Figure Qj] gives some insights into the composi- 
tion of the system at a mixed value of K. At high 
K > 2 the system is completely dominated by the sin- 
gle giant component and there are in fact practically no 
monomers. For N > 2000, K > 1 the number of com- 
ponent clusters is almost identically unity and the num- 
ber of separate monomers almost identically zero. This 
transition remains quite sharp in K even for N > 100 
with a clear gap below the critical K = 1 value. 

We found that the average number of monomers is 
still very small and does not vary with N even for 
K > 1. A fully disconnected system with K = has 
each node as a monomer and the number of monomers 



obviously then grows with N. The first intermediate 
range of < K < Kp shows the system become fully 
connected, and as discussed, Kp — ► 1 as N — ► inf , but 
may be greater than unity for finite N. The next inter- 
mediate range has Kp < K < K c and shows some 
very interesting changes in the system's behaviour. 
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Figure 13: Variation of (modified) mean all pairs dis- 
tance in mixed K networks for N = 75, 100, 140, 
showing maximal characteristic distance at K — 1.5. 



Number of Nodes 

Figure 12: Number of multiple arcs in NK Networks vs 
network size N for K = 1, 2, 3, 4, 5, 6, 7, 8. 



We confirmed empirically our expectation that a sys- 
tem that is constructed to allow self-arcs will exhibit a 
growth in their number linear in N. In the work reported 
in this paper we have typically experimented with sys- 
tems both with and without self- arcs. Figure[l2] shows a 
count of the number of multiple arcs and how these are 
in fact influenced by the presence or absence of self- 
arcs. For the low-K regimes we are most interested in, 
the number of multiple arcs is almost invariant with net- 
work size, although this does grow logarithmically with 
N - particularly with high-K. It appears that neither 
self-arcs nor multiple arcs provide clues to the nature 
of the K = 1.5 circuits behaviour transition. 

Another graph metric that is computationally inex- 
pensive to compute is the all-pairs distance. Elemen- 
tary graph textbooks illustrate this for fully connected 
systems and generally only for bi-directional graphs. 

Figure [13] shows the significance of K = 1.5 and 
needs to be carefully interpreted, recalling that it is mea- 
suring a bulk graph property across regimes where the 
graph is neither necessarily fully connected nor fully 
mutually reachable. As K rises from zero, the N inde- 
pendent nodes become connected and the steepness of 



rise of the all-pairs distance reaches a maximum at the 
percolation transition of K = 1, however even when 
the network forms a single giant cluster, the connec- 
tions are directed and not every node is reachable from 
every other node. Therefore a correct weighted calcu- 
lation of pair-pair distances, where we include two un- 
reachable nodes as contributing zero rather than infinity, 
highlights the structural role of connectivity K = 1.5. 
At this connectivity the network supports the greatest 
set of traversal distances present consistent with being a 
fully connected system. Note further that mutual reach- 
ability tails off rather slowly and even at high K there 
are still "islands of directed disconnection" in the sys- 
tem. 

We can further explore the notion of vertex reacha- 
bility as distinct from connectivity. 

Figure [14] shows the fraction of mutually reachable 
vertices for different mean K in our model, with 560 
nodes. At K < 1 we see multiple disconnected compo- 
nents, but we do not in fact reach full reachability until 
K «6. 

Figure [13] also shows the sharpening of the transition 
with increasing network size. It is worth noting that 
many of the applications of Kauffman networks as they 
relate to real physical and biological systems have very 
definitely finite network sizes N of a few hundred to a 
few tens of thousands and we are not solely interested 
in thermodynamically sized systems. 
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Figure 14: Fraction of Reachable vertices for different 
< K > showing three regimes: multiple components 
for K < 1; one component, but with varying vertex 
reachability for 1 < K < 6 and fully connect and fully 
reachable for large K. 
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Figure 15: Fraction of vertices in the giant component 
for different < K > showing a sharpening of the edge 
for increasing network size N. 



We can verify the bluntening of the percolation tran- 
sition by counting the fraction of the network nodes in 
the giant component. 

Figure [15] shows how our network system genera- 
tion algorithm will produce disconnected clusters even 
above K = Kp = 1 for finite network size N. We con- 
jectured that the K = 1.5 transition might be related 
to partial component disconnection, but even when one 
edits out everything except the giant components and 
arranges them by size, we still see this effect within a 
fully connected component with mean K = 1.5. 



7 Summary and Conclusions 

There seems to be several important transition values 
in the mixed-K model: the percolation transition at 
K w 1 from a disconnected to a fully connected sys- 
tem; the structural transition at K = Ks ~ 1.5 and the 
Kauffman transition from the stable to chaotic regime 
at K c = 2. 

Our study was originally motivated on the assump- 
tion that there is a close relationship between the num- 
ber of elementary circuits in the underpinning structural 
network and the number of attractors that can be sup- 
ported in an associated RBN. As Aldana and Cluzel 
observe [38], the average connectivity i^ e ff appears ir- 



relevant in describing highly heterogeneous scale-free 
topological regimes of an NK system, but it does appear 
to be useful in characterising the transitional regime be- 
tween Kp and K c . 

Specifically, it appears that the structural transition 
at K = 1.5 strongly influences the number of possible 
circuits in the system, but it requires the K = 2 con- 
nectivity transition for the boolean functions to be able 
to exploit it and to produce a number of attractors that 
cross into the chaotic phase. 

In [22] Drossel et al. speculate that the vast number 
of attractors in NK models appear to be a consequence 
of the synchronous updating scheme. While no doubt 
the synchronicity plays a role, we believe our work 
shows that the growth of attractors is also a more fun- 
damental consequence of the structural circuits present 
in the underpinning NK graph. 

We have shown that above Ks the number of struc- 
tural circuits appears to be a relatively simple exponen- 
tial function of the number of connections Na where 
discounting self- arcs and multiple arcs, Na ~ K.N. 
This relationship models both the NK Network system 
and our pair-wise N (K) system. The number of cir- 
cuits Nc therefore does appear to be a useful lower 
bound on the number of attractors in both models, and 
it remains for future work to refine this relationship and 
bound. 
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